Exemplo Vestibular

Author

Luiz Favaro

Curso: MBA DSA USP ESALQ

Fonte: Fávero e Belfiore, MANUAL DE ANÁLISE DE DADOS, Capítulo 09

Prof.: Wilson Tarantin Jr.

Preparação

Instalação e carregamento dos pacotes utilizados

pacotes <- c("plotly", #plataforma gráfica
             "tidyverse", #carregar outros pacotes do R
             "ggrepel", #geoms de texto e rótulo para 'ggplot2' que ajudam a
                        #evitar sobreposição de textos
             "knitr", "kableExtra", #formatação de tabelas
             "reshape2", #função 'melt'
             "misc3d", #gráficos 3D
             "plot3D", #gráficos 3D
             "cluster", #função 'agnes' para elaboração de clusters hierárquicos
             "factoextra", #função 'fviz_dend' para construção de dendrogramas
             "ade4") #função 'ade4' para matriz de distâncias em var. binárias

if(sum(as.numeric(!pacotes %in% installed.packages())) != 0){
  instalador <- pacotes[!pacotes %in% installed.packages()]
  for(i in 1:length(instalador)) {
    install.packages(instalador, dependencies = T)
    break()}
  sapply(pacotes, require, character = T) 
} else {
  sapply(pacotes, require, character = T) 
}
Loading required package: plotly
Loading required package: ggplot2

Attaching package: 'plotly'
The following object is masked from 'package:ggplot2':

    last_plot
The following object is masked from 'package:stats':

    filter
The following object is masked from 'package:graphics':

    layout
Loading required package: tidyverse
Warning in system("timedatectl", intern = TRUE): running command 'timedatectl'
had status 1
── Attaching packages
───────────────────────────────────────
tidyverse 1.3.2 ──
✔ tibble  3.1.8      ✔ dplyr   1.0.10
✔ tidyr   1.2.1      ✔ stringr 1.4.1 
✔ readr   2.1.3      ✔ forcats 0.5.2 
✔ purrr   0.3.4      
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks plotly::filter(), stats::filter()
✖ dplyr::lag()    masks stats::lag()
Loading required package: ggrepel

Loading required package: knitr

Loading required package: kableExtra


Attaching package: 'kableExtra'


The following object is masked from 'package:dplyr':

    group_rows


Loading required package: reshape2


Attaching package: 'reshape2'


The following object is masked from 'package:tidyr':

    smiths


Loading required package: misc3d
Warning: no DISPLAY variable so Tk is not available
Loading required package: plot3D
Loading required package: cluster
Loading required package: factoextra
Welcome! Want to learn more? See two factoextra-related books at https://goo.gl/ve3WBa
Loading required package: ade4
    plotly  tidyverse    ggrepel      knitr kableExtra   reshape2     misc3d 
      TRUE       TRUE       TRUE       TRUE       TRUE       TRUE       TRUE 
    plot3D    cluster factoextra       ade4 
      TRUE       TRUE       TRUE       TRUE 

Se der o erro: “Error in kable_styling”, instale:

sudo apt -y install libfontconfig1-dev e sudo apt install cmake

Análise Exploratória

Carrega a base

load(file = "./data/mod1_01_01_vestibular.RData")

Visualiza dados

Vestibular %>%
  kable() %>%
  kable_styling(bootstrap_options = "striped", 
                full_width = FALSE,
                font_size = 20)
estudante matematica fisica quimica
Gabriela 3.7 2.7 9.1
Luiz Felipe 7.8 8.0 1.5
Patricia 8.9 1.0 2.7
Ovidio 7.0 1.0 9.0
Leonor 3.4 2.0 5.0

Os dados se constituem em 5 alunos e suas médias em 3 matérias: matemática, física e química.

Renomeia linhas

rownames(Vestibular) <- Vestibular$estudante
Warning: Setting row names on a tibble is deprecated.

Visualiza dados em 3D

scatter3D(x=Vestibular$fisica,
          y=Vestibular$matematica,
          z=Vestibular$quimica,
          phi = 0, bty = "g", pch = 20, cex = 2,
          xlab = "Fisica",
          ylab = "Matematica",
          zlab = "Quimica",
          main = "Vestibular",
          clab = "Nota de Quimica")>
  text3D(x=Vestibular$fisica,
         y=Vestibular$matematica,
         z=Vestibular$quimica,
         labels = rownames(Vestibular),
         add = TRUE, cex = 1)

      [,1]  [,2]  [,3]  [,4]
[1,] FALSE FALSE FALSE FALSE
[2,] FALSE FALSE FALSE FALSE
[3,] FALSE FALSE FALSE FALSE
[4,] FALSE FALSE FALSE FALSE

É possível visualizar os clusters formados: a Gabriela e o Ovídio são mais próximos, assim como a Patrícia e o Luiz Felipe. A Leonor está entre os dois grupos.

Estatística descritiva

summary(Vestibular)
  estudante           matematica       fisica        quimica    
 Length:5           Min.   :3.40   Min.   :1.00   Min.   :1.50  
 Class :character   1st Qu.:3.70   1st Qu.:1.00   1st Qu.:2.70  
 Mode  :character   Median :7.00   Median :2.00   Median :5.00  
                    Mean   :6.16   Mean   :2.94   Mean   :5.46  
                    3rd Qu.:7.80   3rd Qu.:2.70   3rd Qu.:9.00  
                    Max.   :8.90   Max.   :8.00   Max.   :9.10  

A média e a mediana das notas de matemática são as maiores, enquanto as da física são as menores. 75% dos alunos tiveram nota de física abaixo de 5, enquanto que em química foram 50% e em matemática 50% tiveram nota a partir de 7.

Padronização

# Se for necessário padronizar, é possível utilizar a função scale()
vest_padronizado <- as.data.frame(scale(Vestibular[,2:4]))
rownames(vest_padronizado) <- Vestibular$estudante
vest_padronizado
            matematica      fisica    quimica
Gabriela    -0.9925328 -0.08223228  1.0369151
Luiz Felipe  0.6616886  1.73373093 -1.1280724
Patricia     1.1055039 -0.66471107 -0.7862322
Ovidio       0.3389136 -0.66471107  1.0084283
Leonor      -1.1135733 -0.32207650 -0.1310387

Com a padronização, todas as matérias ficaram com média zero e desvio padrão igual a 1:

summary(vest_padronizado)
   matematica          fisica            quimica       
 Min.   :-1.1136   Min.   :-0.66471   Min.   :-1.1281  
 1st Qu.:-0.9925   1st Qu.:-0.66471   1st Qu.:-0.7862  
 Median : 0.3389   Median :-0.32208   Median :-0.1310  
 Mean   : 0.0000   Mean   : 0.00000   Mean   : 0.0000  
 3rd Qu.: 0.6617   3rd Qu.:-0.08223   3rd Qu.: 1.0084  
 Max.   : 1.1055   Max.   : 1.73373   Max.   : 1.0369  
sapply(vest_padronizado, sd)
matematica     fisica    quimica 
         1          1          1 

Boxplot

ggplotly(
  Vestibular %>%
    melt() %>%
    ggplot(aes(label = estudante)) +
    geom_boxplot(aes(x = variable, y = value, fill = variable)) +
    geom_point(aes(x = variable, y = value), alpha = 0.5) +
    labs(x = "Variável",
         y = "Nota") +
    scale_fill_manual("Legenda:",
                      values = c("orange", "purple", "bisque4")) +
    theme_bw()
)
Using estudante as id variables
Warning: attributes are not identical across measure variables; they will be
dropped

Aglomeração Hierárquica

Matriz de dissimilaridades (cofenética)

matriz_D <- Vestibular %>% 
  select(matematica, fisica, quimica) %>% 
  dist(method = "euclidean")
matriz_D
          1         2         3         4
2 10.132127                              
3  8.419620  7.186793                    
4  3.713489 10.290287  6.580273          
5  4.170132  8.222530  6.044832  5.473573
data.matrix(matriz_D) %>% 
  kable() %>%
  kable_styling(bootstrap_options = "striped", 
                full_width = FALSE, 
                font_size = 20)
1 2 3 4 5
0.000000 10.132127 8.419620 3.713489 4.170132
10.132127 0.000000 7.186793 10.290287 8.222530
8.419620 7.186793 0.000000 6.580273 6.044832
3.713489 10.290287 6.580273 0.000000 5.473573
4.170132 8.222530 6.044832 5.473573 0.000000

Argumentos para a função dist de cálculo de distância

  • euclidean: distância euclidiana
  • euclidiana quadrática: elevar ao quadrado matriz_D (matriz_D^2)

  • maximum: distância de Chebychev;

  • manhattan: distância de Manhattan (ou distância absoluta ou bloco);

  • canberra: distância de Canberra;

  • minkowski: distância de Minkowski

Clusterização hierárquica

A partir da matriz cofenética, procede-se com a clusterização hierárquica:

cluster_hier <- agnes(x = matriz_D, method = "single")
cluster_hier
Call:    agnes(x = matriz_D, method = "single") 
Agglomerative coefficient:  0.3090455 
Order of objects:
[1] 1 4 5 3 2
Height (summary):
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  3.713   4.056   5.107   5.279   6.330   7.187 

Available components:
[1] "order"  "height" "ac"     "merge"  "diss"   "call"   "method"

Method é o tipo de encadeamento:

  • complete: encadeamento completo (furthest neighbor ou complete linkage)

  • single : encadeamento único (nearest neighbor ou single linkage)

  • average: encadeamento médio (between groups ou average linkage)

Definição do esquema hierárquico de aglomeração

Distâncias para as combinações em cada estágio

coeficientes <- sort(cluster_hier$height, decreasing = FALSE) 
coeficientes
[1] 3.713489 4.170132 6.044832 7.186793

Tabela com o esquema de aglomeração. Interpretação do output:

As linhas são os estágios de aglomeração

Nas colunas Cluster1 e Cluster2, observa-se como ocorreu a junção

  • Quando for número negativo, indica observação isolada

  • Quando for número positivo, indica cluster formado anteriormente (estágio)

  • Coeficientes: as distâncias para as combinações em cada estágio

esquema <- as.data.frame(cbind(cluster_hier$merge, coeficientes))
names(esquema) <- c("Cluster1", "Cluster2", "Coeficientes")
esquema
  Cluster1 Cluster2 Coeficientes
1       -1       -4     3.713489
2        1       -5     4.170132
3        2       -3     6.044832
4        3       -2     7.186793

Números positivos indicam o número do cluster

Números negativos indicam que a observação não entrou em nenhum cluster

Visualização do esquema hierárquico de aglomeração

esquema %>%
  kable(row.names = T) %>%
  kable_styling(bootstrap_options = "striped", 
                full_width = FALSE, 
                font_size = 20)
Cluster1 Cluster2 Coeficientes
1 -1 -4 3.713489
2 1 -5 4.170132
3 2 -3 6.044832
4 3 -2 7.186793

Construção do dendrograma

#dev.off()
fviz_dend(x = cluster_hier)
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> =
"none")` instead.

Dendrograma com visualização dos clusters (definição de 3 clusters)

fviz_dend(x = cluster_hier,
          k = 3,
          k_colors = c("deeppink4", "darkviolet", "deeppink"),
          color_labels_by_k = F,
          rect = T,
          rect_fill = T,
          lwd = 1,
          ggtheme = theme_bw())
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> =
"none")` instead.

Identificando os clusters no dataset

Criando a variável categórica para indicação do cluster no banco de dados

# O argumento 'k' indica a quantidade de clusters
Vestibular$cluster_H <- factor(cutree(tree = cluster_hier, k = 3))
Vestibular
# A tibble: 5 × 5
  estudante   matematica fisica quimica cluster_H
  <chr>            <dbl>  <dbl>   <dbl> <fct>    
1 Gabriela          3.70   2.70    9.10 1        
2 Luiz Felipe       7.80   8       1.5  2        
3 Patricia          8.90   1       2.70 3        
4 Ovidio            7      1       9    1        
5 Leonor            3.40   2       5    1        

Visualização da base de dados com a alocação das observações nos clusters

Vestibular %>%
  kable() %>%
  kable_styling(bootstrap_options = "striped", 
                full_width = FALSE,
                font_size = 20)
estudante matematica fisica quimica cluster_H
Gabriela 3.7 2.7 9.1 1
Luiz Felipe 7.8 8.0 1.5 2
Patricia 8.9 1.0 2.7 3
Ovidio 7.0 1.0 9.0 1
Leonor 3.4 2.0 5.0 1

Estatísticas descritivas dos clusters por variável

ATENÇÃO: Clusters 2 e 3 têm somente uma observação, não calcula sd

  • Estatísticas descritivas da variável matemática

    group_by(Vestibular, cluster_H) %>%
      summarise(
        mean = mean(matematica, na.rm = TRUE),
        sd = sd(matematica, na.rm = TRUE),
        min = min(matematica, na.rm = TRUE),
        max = max(matematica, na.rm = TRUE),
        count = length(matematica))
    # A tibble: 3 × 6
      cluster_H  mean    sd   min   max count
      <fct>     <dbl> <dbl> <dbl> <dbl> <int>
    1 1          4.70  2.00  3.40  7        3
    2 2          7.80 NA     7.80  7.80     1
    3 3          8.90 NA     8.90  8.90     1
  • Estatísticas descritivas da variável ‘fisica’

group_by(Vestibular, cluster_H) %>%
  summarise(
    mean = mean(fisica, na.rm = TRUE),
    sd = sd(fisica, na.rm = TRUE),
    min = min(fisica, na.rm = TRUE),
    max = max(fisica, na.rm = TRUE))
# A tibble: 3 × 5
  cluster_H  mean     sd   min   max
  <fct>     <dbl>  <dbl> <dbl> <dbl>
1 1          1.90  0.854     1  2.70
2 2          8    NA         8  8   
3 3          1    NA         1  1   
  • Estatísticas descritivas da variável química
group_by(Vestibular, cluster_H) %>%
  summarise(
    mean = mean(quimica, na.rm = TRUE),
    sd = sd(quimica, na.rm = TRUE),
    min = min(quimica, na.rm = TRUE),
    max = max(quimica, na.rm = TRUE))
# A tibble: 3 × 5
  cluster_H  mean    sd   min   max
  <fct>     <dbl> <dbl> <dbl> <dbl>
1 1          7.70  2.34  5     9.10
2 2          1.5  NA     1.5   1.5 
3 3          2.70 NA     2.70  2.70

Análise de variância de um fator (ANOVA).

Interpretação do output:

  • Mean Sq do cluster_H: indica a variabilidade entre grupos
  • Mean Sq dos Residuals: indica a variabilidade dentro dos grupos
  • F value: estatística de teste (Sum Sq do cluster_H / Sum Sq dos Residuals)

  • Pr(>F): p-valor da estatística

  • p-valor < 0.05: pelo menos um cluster apresenta média estatisticamente diferente dos demais

  • A variável mais discriminante dos grupos contém maior estatística F (e significativa)

ANOVA da variável ‘matematica’

summary(anova_matematica <- aov(formula = matematica ~ cluster_H,
data = Vestibular))
            Df Sum Sq Mean Sq F value Pr(>F)
cluster_H    2  16.59   8.296   2.079  0.325
Residuals    2   7.98   3.990               

ANOVA da variável ‘fisica’

summary(anova_fisica <- aov(formula = fisica ~ cluster_H,
data = Vestibular))
            Df Sum Sq Mean Sq F value Pr(>F)  
cluster_H    2  32.61   16.31   22.34 0.0429 *
Residuals    2   1.46    0.73                 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

ANOVA da variável ‘quimica’

summary(anova_quimica <- aov(formula = quimica ~ cluster_H,
data = Vestibular))
            Df Sum Sq Mean Sq F value Pr(>F)
cluster_H    2  38.35   19.18   3.506  0.222
Residuals    2  10.94    5.47               

Esquema de aglomeração não hierárquico K-MEANS

Elaboração da clusterização não hieráquica k-means

cluster_kmeans <- kmeans(Vestibular[,2:4],
                         centers = 3)

centers: parametrização da quantidade de clusters

Método de Elbow para identificação do número ótimo de clusters

Apresenta a variação total dentro dos clusters para várias nº de clusters.

Em geral, quando há a dobra é um indício do número ótimo de clusters.

fviz_nbclust(Vestibular[,2:4], kmeans, method = "wss", k.max = 4)

Criando variável categórica para indicação do cluster no banco de dados

Vestibular$cluster_K <- factor(cluster_kmeans$cluster)
Vestibular
# A tibble: 5 × 6
  estudante   matematica fisica quimica cluster_H cluster_K
  <chr>            <dbl>  <dbl>   <dbl> <fct>     <fct>    
1 Gabriela          3.70   2.70    9.10 1         2        
2 Luiz Felipe       7.80   8       1.5  2         1        
3 Patricia          8.90   1       2.70 3         3        
4 Ovidio            7      1       9    1         2        
5 Leonor            3.40   2       5    1         2        

Visualização da base de dados

Vestibular %>%
  kable() %>%
  kable_styling(bootstrap_options = "striped", 
                full_width = FALSE,
                font_size = 20)
estudante matematica fisica quimica cluster_H cluster_K
Gabriela 3.7 2.7 9.1 1 2
Luiz Felipe 7.8 8.0 1.5 2 1
Patricia 8.9 1.0 2.7 3 3
Ovidio 7.0 1.0 9.0 1 2
Leonor 3.4 2.0 5.0 1 2

Análise de variância de um fator (ANOVA)

ANOVA da variável ‘matematica’

summary(anova_matematica <- aov(formula = matematica ~ cluster_K,
                                data = Vestibular))
            Df Sum Sq Mean Sq F value Pr(>F)
cluster_K    2  16.59   8.296   2.079  0.325
Residuals    2   7.98   3.990               

ANOVA da variável ‘fisica’

summary(anova_fisica <- aov(formula = fisica ~ cluster_K,
                            data = Vestibular))
            Df Sum Sq Mean Sq F value Pr(>F)  
cluster_K    2  32.61   16.31   22.34 0.0429 *
Residuals    2   1.46    0.73                 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

ANOVA da variável ‘quimica’

summary(anova_quimica <- aov(formula = quimica ~ cluster_K,
                             data = Vestibular))
            Df Sum Sq Mean Sq F value Pr(>F)
cluster_K    2  38.35   19.18   3.506  0.222
Residuals    2  10.94    5.47               

Comparando os resultados dos esquemas hierárquico e não hierárquico

Vestibular %>%
  select(estudante, cluster_H, cluster_K) %>%
  kable() %>%
  kable_styling(bootstrap_options = "striped",
                full_width = FALSE,
                font_size = 20)
estudante cluster_H cluster_K
Gabriela 1 2
Luiz Felipe 2 1
Patricia 3 3
Ovidio 1 2
Leonor 1 2

FIM!